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Abstract 

The present contribution proposes a general constitutive model to simulate the or¬ 
thotropic stiffness, pre-peak nonlinearity, failure envelopes, and the post-peak softening 
and fracture of textile composites. 

Following the microplane model framework, the constitutive laws are formulated 
in terms of stress and strain vectors acting on planes of several orientations within 
the material meso-structure. The model exploits the spectral decomposition of the 
orthotropic stiffness tensor to define orthogonal strain modes at the microplane level. 
These are associated to the various constituents at the mesoscale and to the material 
response to different types of deformation. Strain-dependent constitutive equations are 
used to relate the microplane eigenstresses and eigenstrains while a variational principle 
is applied to relate the microplane stresses at the mesoscale to the continuum tensor at 
the macroscale. 

Thanks to these features, the resulting spectral stiffness microplane formulation can 
easily capture various physical inelastic phenomena typical of fiber and textile compos¬ 
ites such as: matrix microcracking, micro-delamination, crack bridging, pullout, and 
debonding. The application of the model to a twill 2x2 shows that it can realisti¬ 
cally predict its uniaxial as well as multi-axial behavior. Furthermore, the model shows 
excellent agreement with experiments on the axial crushing of composite tubes, this 
capability making it a valuable design tool for crashworthiness applications. 

The formulation is computationally efficient, easy to calibrate and adaptable to 
other kinds of composite architectures of great current interest such as 2D and 3D 
braids or 3D woven textiles. 


1. Introduction 

Thanks to their excellent specific mechanical performances and the recent develop¬ 
ments in manufacturing technologies, the range of engineering applications of textile 
composites is continuously expanding. Current applications include land, marine and 
air transportation, wind and tidal energy production, and blast protection of civil in¬ 
frastructures and vehicles nm. However, in order to take advantage of the outstand¬ 
ing characteristics of these materials, design tools to simulate the orthotropic stiffness, 
pre-peak nonlinearity, failure envelopes, and the post-peak softening and fracture are 
quintessential. 
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Since the pioneering works by Ishikawa and Chou m E] and Ishikawa et al. [B], 
several formulations have been proposed, with varying degrees of success, to model 
the elastic properties of textile composites [M3] and their failure mechanisms m- 
[20] • In general, however, these models stand on strength criteria to describe failure of 
the mesoscale constituents thus lacking completely of any description of the fracture 
mechanics involved. This is a serious dehciency being extensive intra-laminar cracking 
one of the main failure mechanisms in most applications of textile composites. 

Modeling the fracturing behavior of textile composites, not only requires a fracture 
mechanics framework, it also urges the acknowledgment of their quasi-brittle character 
which highly affects the process of crack nucleation and growth. In facts, due to the 
complex mesostructure characterizing quasi-brittle materials (such as composites and 
nanocomposites, ceramics, rocks, sea ice, bio-materials and concrete, just to mention a 
few), the extent of the non-linear Fracture Process Zone (FPZ) occurring in the pres¬ 
ence of a macrocrack is usually not negligible [21]. The stress held along the FPZ is 
nonuniform and decreases with crack opening gradually, due to discontinuous cracking, 
crack bridging by hbers, and frictional pullout of inhomogeneities. As a consequence, 
the fracturing behavior and, most importantly, the energetic size effect and the qua¬ 
sibrittleness effect associated with structure geometry, cannot be described by means 
of the classical Linear Elastic Fracture Mechanics (LEFM). To capture the effects of a 
hnite FPZ size, the introduction in the formulation of a characteristic (hnite) length 
scale of the material is necessary [211122] . This is attempted in the present work. 

Inspired by a recent theoretical framework for unidirectional composites by Cusatis 
et al. [23] El] , this contribution aims at proposing a general constitutive model to simu¬ 
late the damaging and fracturing behavior of textile composites. The formulation stands 
on the dehnition of strain-dependent constitutive laws in terms of stress and strain vec¬ 
tors acting on planes of several orientations within the material meso-structure. In 
this way, the model can easily capture various physical inelastic phenomena typical of 
hber and textile composites such as: matrix micro cracking, micro-delamination, crack 
bridging, pullout, and debonding. 
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Thanks to the coupling with the crack band model |2Sl I2S], the formulation is 
endowed with a characteristic length dependent on the strength and the fracture energy 
of the material. This is key to capture the intra-laminar size effect, a salient feature of 
composite structures. This aspect, too often overlooked in the literature on composites, 
is a determinant factor for damage tolerance design of large composite structures. 

2. Theoretical framework 

2.1. Microplane model 

Inspired by the slip theory of plasticity pioneered by Taylor 1271 and later rehned 
by Batdorf and Budiansky [28], the microplane theory was originally developed to 
describe the softening damage of heterogeneous but statistically isotropic materials 
such as concrete and rocks [291 EQ|. Since its introduction in the early 1980s, the 
microplane model for concrete has evolved through 7 progressively improved versions 
labeled as Ml [221130], M2 [21], M3 [22], M4 [221 [21], MS [2^, M6 [22], M7 [27] and it 
has been recently adopted for the simulation of concrete at early age [2E]. Microplane 
models have also been developed for other complex materials such as jointed rock [29] . 
sand, clay, rigid foam, shape memory alloys, and unidirectional and textile composites 
[221122] 1201122] . A high order microplane model [22] was also derived recently on the 
basis of an underlying discrete model [IG] [27] . 

A key feature of the microplane model is that the constitutive laws are formulated in 
terms of the stress and strain vectors acting on a generic plane of any orientation within 
the material meso-structure, called the microplane. These planes can be conceived as 
the tangent planes of a unit sphere surrounding every point in the three-dimensional 
space (Fig. [^). The microplane strain vectors are the projections of the macroscropic 
strain tensor, whereas the macroscopic stress tensor is related to the microplane stress 
vectors via the principle of virtual work. The adoption of vectors rather than tensors 
makes the approach conceptually clearer while the introduction of microplanes allows 
to inherently embed the effect of the mesostructure into the formulation. 
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In this contribution, a kinematically constrained microplane model is adopted. This 
means that the strain vector on each microplane is the projection of the macroscopic 
strain tensor. In kelvin notation |1HI SH] this reads: 


£p = Pe 


( 1 ) 


where Sp = £l\^ represents the microplane strain vector (Fig [^) with = 

normal strain component and Em and Sp = shear strain components. Further, 
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is a 3 X 6 matrix relating the macroscopic strain tensor to the microplane strain as a 
function of the plane orientation. In facts, Nij = UiUj, Mij = {miUj + mjni)/2 and 
Lij = {liUj + ljni)/2, where rij, and U are local Cartesian coordinate vectors on the 
generic microplane with rij being the i-th component of the normal (Fig. With 

reference to the spherical coordinate system represented in Fig. &). the foregoing 
components can be expressed as a function of the spherical angles 9 and (p: ni = 
sin 0 cos (p, 712 = sin 0 sin (p, = cos 6 while one can choose mi = cos 0 cos (p, m 2 = 

cos 0 sin (p, m 3 = — sin 0 which gives, for orthogonality, li = — sinyj, I 2 = cos ip and 

h = 0 . 

According to the microplane framework, the constitutive laws are then dehned at 
the microplane level in a vectorial form. This makes the formulation conceptually clear 
and allows embedding the effect of the direction of damage in the constitutive law 
automatically. After the microplane stress vectors crp are computed, the macroscopic 
stress tensor is dehned in a variational sense through the principle of virtual work: 


<T = - / P’^CTpdVt 

27r ./o 


(3) 


where is the surface of a unit sphere representing all the possible microplane orien¬ 
tations. 
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2.2. Spectral decomposition of the elastic tensor 

In the microplane formulation, the material anisotropy is addressed by decomposing 
the stress and strain tensors into energetically orthogonal modes through the spectral 
stiffness decomposition theorem [SOHSS]. The following sections are intended to provide 
a brief introduction of the theory. 

2.2.1. Spectral decomposition of the elastic tensor 

The elastic behavior of a general anisotropic material can be expressed in Kelvin 
notation [IHl 09] as: 

a = Ce (4) 

where cr = [an a 22 cr^s \/2a23 \/ 2 cti 3 V^a^] , £ = [^n £22 ^33 ^^£23 
are the contracted forms of the stress and strain second-order tensors and C represents 
the contracted form of the fourth-order elastic tensor. The indices refer to Cartesian 
coordinates Xi {i = 1,2,3) as dehned in Figures ([^,b). It is worth mentioning here 
that the factor \/2 assures that both the stiffness tensor and its column matrix have 
the same norm, given by the sum of the squares of their elements. 

According to the spectral decomposition theorem [50II53] . the stiffness matrix C 
can be decomposed as follows: 

c = Y. ( 5 ) 

/ 

where / = 1, 2...6, are the eigenvalues of the stiffness matrix and 

are a set of second-order tensors constructed from the elastic eigenvectors <f>j. The J-th 

eigenvector cpj has multiplicity n and is normalized such that 

In general, the elastic eigenvalues and eigenmatrices can be found by solving the 
following eigenvalue problem: 

[C - v = 0 (6) 

where, it should be noted, the number of independent elastic eigenvalues is strictly 
related the degree of anisotropy of the material. In general anisotropy, the solution of 
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Eq. (|^ provides 6 independent eigenvalues which are functions of the 21 independent 
elastic constants of the material. 

It is worth mentioning that the eigenmatrices are: 

1 . partition of unity, i.e. 

2. orthogonal, i.e. = 0 if / 7 ^ J; 

3. idempotent, i.e. 

Moreover, it should be remarked that the property 2) above leads to energetic orthog¬ 
onality: dWij = a-JSsj = = 0 . 

An important characteristic of the elastic eigenmatrices is that they provide 
a way to decompose the stress and strain tensors into energetically orthogonal modes. 
These are called here eigenstresses and eigenstrains and are dehned as: 

(Tj = and £/ = (7) 

It is easy to show that cr = 'Yli ^ind e = 'Yli whereas the relation between 
eigenstresses and eigenstrains can be found introducing the related elastic eigenvalues: 

CTi = 

A similar analysis can be applied to the compliance elastic matrix S = C~^. In this 
case, the relation between the compliance matrix and its eigenmatrices reads: 

5 = ( 8 ) 

1 

where it could be easily shown that the compliance eigemenmatrices coincide with the 
ones derived from the stiffness tensor, 

2.2.2. Microplane model with spectral decomposition 

In this section, the decomposition of the strain tensor into energetically orthogo¬ 
nal eigenstrains through the spectral decomposition of the elastic tensor is adopted 
to extend the microplane approach to general anisotropic materials. By the spectral 
decomposition of the strain tensor and a separate projection of each eigenstrain, each 
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microplane vector can be decomposed into microplane eigenstrain vectors as: 


N 

Sp = where 

/ 


eP = PeU) = 
p(b = pc^i) 


( 9 ) 


where N = number of independent eigenmodes. 

The main advantage is that now, different constitutive laws describing the material 
behavior at the microplane level can be related to each eigenmode. As it will be clear 
in the following sections, this not only allows the description of the material anisotropy 
but also to address the different damaging mechanisms related to different loading con¬ 
ditions. Accordingly, from the microplane eigenstrains, the microplane eigenstresses 
crP can be dehned through specihc constitutive laws: crP = f (spi, Sp 2 ...) sP. Then, 
substituting the microplane eigestresses into Eq.([^ and recalling Eq. Q, the macro¬ 
scopic stress tensor can be computed in a variational sense through the principle of 
virtual work [23]: 


N 


<T = 


27r 


'O r 


crPdQ 


( 10 ) 


where is the surface of a unit sphere representing all the possible microplane ori¬ 


entations. It is worth mentioning here that Eq. (10) represents a weak variational 
constraint. In general, the projection of the stress tensor does not coincide with the 
microplane eigenstress, i.e. crP ^ P^^'^cr. Such a coincidence, called double constraint, 
holds in the elastic regime if and only if the microplane eigenstress vector is proportional 
to microplane eigenstrain vector through the related eigenvalue, as shown by Cusatis 
et al. 


aP = X^eP. 


2.3. Analysis of microplane eigenstrain modes and physical interpretation 

It is useful to study the distribution of the normal components of each microplane 
eigenstrain on the microplane sphere to give a physical interpretation. This is key 
to dehne physically-based constitutive laws capable of addressing both material and 
damage-induced anisotropy. 
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2.3.1. Strain decomposition for a twill 2x2 composite 

In order to exemplify the application of the spectral stiffness decomposition theorem 
to the microplane model, let us consider the case of a textile composite reinforced by 
a 2x2 twill fabric made of carbon hbers. The elastic properties of the material can 
be found in Table Q while the Representative Unit Cell (RUC) for this material is 
represented in Figure ([^) together with the local coordinate system used to define the 
elastic tensor of the material. As can be noted, planes 1 - 2 , 1-3 and 2-3 represent 3 
planes of material symmetry so that the material can be considered as orthotropic [M] • 
For an orthotropic material, the number of independent elastic constants is 9. How¬ 
ever, for the case of a twill 2 x 2 , the in-plane elastic moduli as well as the out-of-plane 
Poisson ratios and shear moduli are usually almost the same (i.e. k. E 2 = E, 
z /13 z /12 = z/' Gi 3 ~ Gi 2 = G' respectively). Accordingly, the elastic compliance 

tensor. S', can be expressed as a function of only 6 elastic constants as: 
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where u = z /23 = 1^32 = in-plane Poisson ratio and G = G 23 = in-plane shear modulus. 
The elastic stiffness tensor can be computed as C = S~^. 

Thanks to the material symmetry, the solution of the elastic eigenvalue problem 
in Eq. (|^ can be carried out analytically giving the following 5 independent elastic 




eigenvalues: 
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( 12 ) 


The related eigematrices, constructed from the normalized eigenvectors, are: 

,,(1) _ X^^lj + X + <^ljF + + ^2j {^12 + ^Si) + ^3j (<^i3 + ^ 21 ) 


2G' 


c)r = 




2 + x^ 


^(2) _ ((^2j - hj){5i2 - 5iz) 

ij 2 


(3) _ ^ {SljS2i + ^2j^lt + Sijd^i + dsjdii) + 62 j {5i2 + ^3i) + ^3i {Si3 + S2i) (13) 


cr = 




2 +e 
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C\f = <5,454i 




where 5ij = Kronecker delta, f,j = 1,2 ...6 and: 


X = 
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(14) 


are dimensionless constants which depend on the elastic properties of the material 
system. 

2.3.2. Microplane eigenstrain distribution 

After the spectral decomposition of the stiffness matrix, the expressions for the nor¬ 
mal, and total shear eip'^ = \J{e^mY + components can be easily computed 

substituting Eqs.([^ and (13r-e) into Eq. (|^. As can be noted from their expres¬ 
sions, reported in Table (|^, the microplane eigenvector components depend on 1) the 
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the applied macroscopic strains through the functions a = (£2 + £3 + ^ix) /(2 + 
j3 = 1/2 (£3 — £2), 7 = (^2 + ^3 + £iC) /(2 + 2) the orientation of the microplane 

normals and 3) the material elastic properties through the dimensionless constants ^ 


and X dehned in Eqs.(14i,b). 


According to Eqs. (14i,b), ^ = —10.0 and x = 0.2 respectively for the material 
under study and the distribution of normal and shear components can now be analyzed 
for different applied strains at the macroscale. Let us now consider the distribution of 
the normal strain components on the microplane sphere caused by a macroscale uniaxial 
strain applied along the X2-axis or the xa-axis (Fig. & ,b). Strain mode 4 and strain 
mode 5 are exactly zero because they do not depend on 82 and £3 whereas mode 3 is 
nonzero but negligible compared to modes 1 and 2. As can be noted from Figure ([^, 
which represents the distribution of the normal components of mode 1 and 2 on the 
microplane sphere, mode 1 acts as a volumetric-like mode, being uniformly distributed 
along the 2-3 plane. Conversely, mode 2 acts as a deviatoric-like mode, loading mainly 
the microplanes whose normal is close to the direction of application of the macroscale 
strain. The two modes have almost the same magnitude in the direction of application 
of the macroscale strain (Figure |^,b). Since both modes describe the uniaxial behavior 
of the material system in the plane of the fabric, it is convenient to consider the effect 
of the two modes summed together: = \i/\ 2 e^p + e^p\ The resulting mode will 

be referred as mode 12 in the following. It should be noted that mode 1 is weighted by 
the ratio between the hrst and the second eigenvalues of the elastic tensor, which for 
most of textile composites is within 0.8 — 1.2. The reason of this choice will be clear 
in section 2.4.3 As Figure shows, mode 12 describes the in-plane uniaxial behavior 
of the material in the direction of application of the load. Being the in-plane behavior 
a “hber-dominated” property, this analysis shows that the response of the composite 
material subjected to strain mode 12 strongly depends on the behavior of the hbers. 

A different picture arises if one considers a uniaxial strain at the macroscale applied 
in the out-of-plane direction, i.e. xi-axis. Strain mode 4 and strain mode 5 are exactly 
zero because they do not depend on Si while mode 12 is nonzero but negligible compared 
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to mode 3. Figure shows the distribution of the normal strain component. It can be 
concluded that mode 3 is mainly related to the behavior of the material system when 
subjected to an out-of-plane uniaxial strain at the macroscale, ei. The behavior in this 
direction is expected to be “matrix-dominated” as the fabric and the matrix mainly act 
in series coupling in this direction of loading. 

As can be expected from the equations describing its microplane strain components 
(Table [^, mode 4 is activated only in the presence of a macroscopic in-plane shear 
strain, i.e. £ 4 . In this conhguration, all the other modes are zero while mode 4 acts 
mainly on microplanes oriented in direction 6 = 7 r/ 4 . The in-plane shear behavior is 
highly dictated by the matrix. Accordingly, as for mode 3, mode 4 can be considered 
as a “matrix dominated mode”. Similar considerations can be drawn for mode 5, this 
mode being activated in the presence of macroscopic out-of-plane shear strains, i.e. £5 
and £ 0 . For this case, the distribution of the normal component of the microplane 
eigenstrain is shown in Figure]^ 

2.4- Formulation of microplane constitutive laws 

In this section, the microplane constitutive laws which, for each eigenmode, provide 
the relationship between the strain vector and the stress vector at the microplane level, 
are presented. 


2.4-1- Elastic behavior 

The elastic behavior is formulated by assuming that normal and shear eigenstresses 
on the microplanes are proportional to the corresponding eigenstrains: 


_ ^{I) {!) (I) _ x(/) 

O N — £jy, O — A £ 


(T 

M ) 


a 


(I) 


= 


-T) 


(15) 


where = J-th elastic eigenvalue. It should be noted that, as proved by Cusatis et 
al. in [23], these relations guarantee a double constraint in the elastic regime, i.e. the 
projection of the stress tensor does coincide with the microplane eigenstress: cTp^ = 
P^^'^cr. This condition is not generally true and it is violated as soon as the material 
reaches the inelastic regime. 
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2.4-2. Inelastic behavior: general constitutive law at microplane level 

The spectral stiffness microplane formulation of the nonlinear and inelastic behavior 
aims at representing the main physical mechanisms characterizing the meso-scale failure 
such as e.g. fiber failure, fiber-matrix debonding and matrix microcracking. 

Let us consider a generic J-th eigenmode and, following Cusatis et al. [16], [55] , dehne 
a microplane effective eigenstrain as: 




( 16 ) 


where ^ total shear strain component of J-th microplane eigen¬ 

strain. The constitutive law can then be dehned by means of an effective eigenstress, 
The relation between the stress and strain microplane components can be found 
imposing the consistency of the virtual work: 

r(T 


SWj — {e]\[6e]\[ -IsmSsm + + + 

(17) 


where, it should be noted, thanks to the dehnition in Eq. (16) and to Eq. dl?] ), dVdi > 0 
for all the microplanes. This means that energy dissipation on each microplane is non¬ 
negative, a sufficient condition to fulhll the second law of thermodynamics. 

By means of Eq. the relationship between normal and shear stresses versus 

normal and shear strains can be formulated through damage-type constitutive equa¬ 
tions: 

= = ( 18 ) 


= [a 


It is worth pointing out that, substituting Eqs. (18) into ([16|) and rearranging, leads 


to the following dehnition of the J-th effective stress: where 

(jy^ = \J{cr^My -|- {o'^L^y= total shear stress component. 

The effective stress is assumed to be incrementally elastic, i.e. and 

it is formulated such that 0 < where 

with subscript i =t for tension and i =c for compression is a limiting boundary enforced 
through a vertical (at constant strain) return algorithm. It is worth mentioning here 
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that, in general, is a fnnction of the microplane orientation and of the eqnivalent 
strains pertaining to other modes. This allows to inherently embed in the formnlation 
the effects of damage anisotropy and the interaction between damaging mechanisms. 


As clarihed in section 2.4.7 each mode can be ascribed to a particnlar micromechan¬ 
ical constitnent and a type of deformation. This make them particnlary snitable for 
the description of the diverse damage mechanisms of the material throngh eigenmode- 
dedicated bonndaries, the description of which being the snbject of the following sec¬ 
tions. 


3.4-3. Inelastic behavior: mode 12 


As discnssed in section 2.3, since both modes 1 and 2 describe the hber behavior 


in the presence of a macroscopic nniaxial in-plane strain, it is convenient to combine 
them together. Accordingly, the formnlation of the effective strain and stress has to be 
changed to gnarantee work consistency. 

Let ns dehne mode 12 effective microplane eigenstrain as: 


^( 12 ) ^ 


A(b r 


A( 2 ) 






( 2 )^ 


-( 2 )^ 


(19) 


After introdncing the eqnivalent stress, crd^)^ the virtnal work can be calculated as done 


in Eq. (17): 


5>Vi2 = 
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( 20 ) 

Now, recalling that (5>Vi2 = X)/=i relationship between 

normal and shear stresses versus normal and shear strains can be formulated through 
the following damage-type constitutive equations: 
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£( 12 ) 
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L 


( 21 ) 


Again, combining Eqs. (|21|) with 
a(i 2 ) = 


one gets the new dehnition of the effective stress: 


)^^) = \J A(2)/A(i)[((j^^)2 -I- (cr^^^)2] -I- (crj^^)2 -I- (cr^^^)2. The aforementioned dehnition 
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of the effective strain and stress guarantees the correct description of the elastic behav¬ 
ior. As a matter of fact, assuming the effective stress to be incrementally elastic, 
i.e. one obtains again Eqs. (15). 


The strain dependent limiting boundary, ...^O^ip), is expressed as 

follows: 
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( 22 ) 


where the brackets (•) are used in Macaulay sense: (x) = max(a;, 0). The functions 
s^^'^\6,ip,£^mL) = and c^^‘^\9,ip,£ml^) = represent mode 12 mi¬ 

croplane tensile and compressive strength respectively which, in general, depend on the 
microplane orientation. However, for the material under study and for most of sym¬ 
metric textile composites in which the behavior in the warp direction is similar to the 
weft direction, they can be both assumed constant. 

The microplane strengths depend, in general, on the maximum effective mode 4 
strain, £mlx. The interaction between the different modes will be clarihed in section 

Mm 

As can be seen in Figure ([^), the boundary evolves exponentially as a func¬ 
tion of the maximum effective strain, which is a history-dependent variable dehned as 
= max.r<t[£^^^^^(^)] with i =t for tension and i =c for compression. The ex¬ 
ponential decay of the boundary alf starts when the maximum effective strain reaches 
its elastic limit el^^^\9,ip,emlx) = and el^^\9,ip,emlx) = ten¬ 

sion and compression respectively. The decay rate is governed by the post-peak slope 
(Figure]^) which, as can be derived from Eq. (22), reads: 


fT (A12)\ _ -,^^<^.12 

Xlt \ fcmax I — (12) 1 (12) I 

^bt \ ^bt J 

-( 12 ) ^ f 

■ C> I ) 


atl2 —1 


in tension 


acl2 —1 


(23) 


Hr I £max I = "bir I — .5,1 ' — I in compression 

For ati 2 ,aci 2 < 1 the initial softening modulus goes to inhnity, i.e. the curve has an 
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initial vertical slope. This seems to match perfectly the in-plane behavior of composites 
which behave elastically up to the peak in stress before softening. In general, the total 
number of required parameters to the describe mode 12 in tension and compression is 
6. However, in absence of specihc experimental data, the exponents ati 2 , aci 2 can be set 
equal, as done in the present work (see Table [^. 


2.4-4- Inelastic behavior: mode 3 

As shown in the previous sections, the 3rd mode can be related to the normal defor¬ 
mation in the out-of-plane direction. The strain dependent boundary can be expressed 
by the following equations: 


ajf = s(3) (p) for sm > 0 

a) for £jv3 < 0 


(24) 


where { 9 , ip) = = mode 3 microplane tensile strength, ( 9 , (p) = = mode 3 

microplane compressive strength. These functions are both assumed to be independent 
of the microplane orientation. As can be seen in Figure (|^), Eqs. (24i,b) result into a 
bilinear boundary in tension and compression. The parameters and kic^ represent 
the value of the effective strain at the beginning of the linear softening region whereas 
the parameters k^f and kl^^ are defined such that k^^^ + kj^f and k^^^ -h kic^ represent the 
value of the effective strain when = 0 in tension and compression respectively. The 
total number of required parameters to describe mode 3 in tension and compression is 6. 
However, in most of applications of composite laminates, the out-of-plane stress is often 
negligible. Accordingly, the material is likely to behave elastically in this direction so 
that an accurate characterization of the parameters related to mode 3 is not required. 
In the present work, in the absence of direct experimental data, mode 3 parameters are 
estimated from [SUll37j assuming that the out-of-plane behavior is “matrix-dominated”. 


2.4-3. Inelastic behavior: mode 4 

The 4th mode is related to in-plane shear deformation. In this type of deformation, 
the material is subjected to extensive microcracking before reaching the peak in stress 
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which produces a remarkable non-linear behavior. In order to capture this phenomenon, 
the strain dependent boundary is expressed by the following equations: 
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for < 0. The functions si‘*i {0, ^p) — s[f^ and ci‘** {0, p) — represent mode 
4 microplane tensile and compressive stresses at which non-linearity starts and they 
are both assumed to be independent of the microplane orientation, p is a parameter 
controlling the slope of the curve in the pre-peak region and it assumed to be the same 
in tension and compression. The parameters k^ac represent the value of the 

effective strain at the beginning of the linear softening region whereas the parameters 
k^ and k^^^ are defined such that k^ + k^^t k^^^ -|- k^ac represent the value of 
the effective strain when = 0 in tension and compression respectively. A typical 
boundary resulting from Eqs. (25) and (26) can be seen in Figure (|^). The total 
number of required parameters to describe mode 4 in tension and compression is 7. In 
cases in which the macroscopic in-plane shear behavior is independent of the sign of 
the deformation (which is not always true for textile composites), the same parameters 
at microplane level can be used in tension and compression. In such a case, considered 
in the present work, the number of required parameters is 4 as reported in Table 


2.4-6- Inelastic behavior: mode 5 

Mode 5 describes the shear deformation in the out-of-plane direction. In the ab¬ 
sence of direct experimental data, it is assumed here that the constitutive behavior 
at microplane level resembles the one related to mode 4. It is worth noting that the 
effects of these assumption are considered to be negligible since, in practical situations. 
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out-of-plane shear deformation will always been within the elastic regime. In facts, 
inter-laminar cracks would emerge far before than reaching any nonlinear behavior in 
the material. 


Interaction between eigenmodes 

After dehning the inelastic boundaries, the interaction between modes needs to be 
clarihed to capture the complex multi-axial behavior of textile composites. 

The dehnition of mode 12 effective strain accounts already for the interaction be¬ 
tween the volumetric and the deviatoric in-plane components at microplane level. As 
a matter of facts, the stress peak for mode 12 is reached when with the 

subscript i = t for tension and i = c for compression. This corresponds to the following 
microplane strain multi-axial criterion: 
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where the dehnition of effective strain reported in Eq. (15) has been used for mode 1 


and 2. As soon as Eq. (27) is satished, strain softening starts to occur. 

The interaction between mode 12 and mode 4 can be described in a similar way: 

(28) 
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where again i = t for tension and i = c for compression. As soon as Eq. (28) is satished. 


the new strains at peak are updated and strain softening starts to occur. This corre¬ 
sponds to a vertical scaling of the softening boundary of a factor of \/l 
for mode 12 and 
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values of the ehective strain for mode 4 and 12 respectively as soon as Eq. (28) is 
fulhlled. 

Considering the deformations associated to mode 12 and mode 4, it is easy to 


understand that the criterion in Eq. (28) dehnes the mechanical behavior of the material 
subjected to simultaneous normal and shear in-plane deformations. 
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Similar interaction laws for the other modes could also be defined. However, in the 
absence of experimental data for validation, only the interaction between modes 1, 2 
and 4 has been considered in the present work. 

2.4-8. General unloading-reloading rule at microplane level 

The inelastic formulation needs to be completed by unloading-reloading rules to 
simulate cycling loadings and unloading of the material next to a damaged zone un¬ 
dergoing softening. Figure ([^) illustrates the unloading-reloading rule adopted in this 
work in terms of effective stress versus effective strain. 

Let us assume that unloading occurs after the effective strain increased continu¬ 
ously from zero to a certain value, eiiax where index i =t for tension and i =c for 
compression. The effective stress decreases elastically until it reaches a zero value and 
remains constant at zero for further decreases of the effective strain. During reloading, 
the effective stress remains zero until the effective strain reaches the reloading strain 
limit, and, beyond this point, the behavior is incrementally elastic. The reloading 

strain limit is dehned as where is assumed to be a 

material parameter which governs the size of hysteresis cycles and, consequently, the 
amount of energy that the material can dissipate during cycling loading. For = 1, 
the dissipated energy density is zero, whereas for = 0, it is a maximum and equal 
to alPsm^x ■ In the present contribution, in the absence of specihc experimental data, 
the dissipated energy density was set equal to zero for all the eigenmodes. 

3. Materials and tests for calibration and validation 

3.1. Materials 

In order to calibrate and validate the model, experiments were conducted on twill 
2x2 woven composite specimens manufactured by compression molding. A DGEBA- 
based epoxy resin was chosen as polymer matrix whereas the reinforcement was provided 
by a twill 2x2 fabric made of carbon hbers. A constant thickness of approximately 1.9 
mm, corresponding to 8 laminae, was used for all the tests. Two different lay-ups. 
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namely a [0°]8 and a Quasi-Isotropic (Q.I.) [0°/45°/—45°/90°]s, were tested in order to 
provide sufficient data for calibration and validation. 

3.2. Uniaxial tests 

Uniaxial tests on coupons were conducted to obtain the elastic properties and 
strength of the composite. In order to provide sufficient data in tension and com¬ 
pression and to provide information regarding the shear behavior of the material, the 
following test configurations were adopted: 

1. 0° conhguration (conducted on both lay-ups): where the loading axis is aligned 
with direction 3 of the RUG (see Fig. [^, to yield axial moduli and strength of 
the composite lamina, both under tension and compression; 

2. 45° configuration (conducted for only the [0°]8, under tension), where the loading 
direction is at an angle of 45° with direction 3 (see Figure [^, to yield the shear 
modulus and strength of the lamina. 

The uni-axial tension tests on the [0°]8 lay-up were used for calibration of the model 
whereas, those on the QI layup were used for validation (i.e. no parameter of the 
model was changed to match the experimental results). All the tests were conducted 
in accordance with ASTM standards [SSI ES] • 

3.3. Size effect tests 

Following Bazant et al. [SU [22], size effect tests were conducted on single-edge- 
notched tension (SENT) specimens, using the [0°]8 lay-up. These tests represent an 
indirect way to characterize the initial inta-laminar fracture energy of the material, 
Gf, which is usually hard to measure directly through uniaxial tests due to snap-back 
instability [221 |2Sl SSI ED] • 

The tested SENT coupons were of three different sizes, scaled as 1:2:4. The dimen¬ 
sions of the smallest coupon were L =100 mm, H = 20 mm and T =1.9 mm. The length 
L and width D of the medium and large specimen were scaled accordingly while the 
thickness T remained constant. Each coupon had on one edge a pre-cut notch of length 
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a = 0.2D. More details on these tests will be provided in a forthcoming contribution 

ra. 

3.4- Impact tests on composite crash can 

Finally, to verify the predictive capability of the model, impact tests on composite 
crash cans were conducted. The crush can consisted of two parts manufactured by com¬ 
pression moulding, namely a hat section tube and a reinforcing plate glued together by 
toughened epoxy glue along flanges (Figure [T^,d). Two different lay-up configurations 
were studied, namely: 

1. hat section tube: [0°]ii, plate: [0°]8; 

2. hat section tube and plate: [07907457-4570790707-4574579070°]; 

The composite tubes, accurately fixed at the bottom, were impacted by a flat mass of 
74.4 kg at a velocity of 4.6 m/s in a drop tower. 

4. Calibration and validation 
4.I. Finite element implementation 

The present model was implemented in Abaqus Explicit v6.11 as a user material 
subroutine VUMAT [62] and in the special purpose software MARS [63|. The struc¬ 
tures studied were meshed using four-node shell elements with reduced integration and 
hourglass control [M] ES] • 

4-1.1. Crack band model for strain localization 

To ensure objective numerical results in the presence of strain localization, the crack 
band model proposed by Bazant et al. [25l |26l |66| is adopted. In this approach, the 
width of the damage localization band, Wc, is considered as a material property. This 
width is also equal to the mesh size, which was here chosen as 2 mm. A change in the 
element size requires the scaling of the post-peak response of the material such that the 
fracture energy remains unchanged. The size Wc should not be confused with the size 
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of the RUC, which is roughly 7.5 mm and merely represents the repeating geometric 
unit of the material. 

Thanks to the crack band model, a characteristic size of the material is inherently 
embedded in the formulation. This is key to correctly describe the intra-laminar frac¬ 
turing of the material as will be clear in the following sections. 


4-2. Strength and post-peak fracturing behavior 

For the calibration of the model, it is now possible to take advantage of the physical 


interpretation given to each eigenmode in section |2.3[ It is known that the strength 
and fracture behavior in the 0° conhguration is governed by mainly the properties of 
the fabric, and, in the 45° conhguration, mainly by the matrix properties in shear. 

It was shown that mode 12 is the “hber-dominated” mode governing the in-plane be¬ 
havior under normal loading condition whereas mode 4 is the “matrix-dominated” mode 
describing the shear behavior. Accordingly, being these modes completely independent, 
the tests in 0° conhguration can be used to uniquely characterize the constitutive law 
related to mode 12 whereas the tests in 45° conhguration can be considered for mode 
4. Thus, thanks to the clear physical interpretation of each mode, the calibration of 
constitutive laws becomes straightforward. 

In general, the calibration is performed in two steps: (1) Determination of parame¬ 
ters from computations for one material point and (2) verihcation of structural response 
from coupon level simulations. The need of a two-step process is explained as follows. 
As far as the elastic stress-strain behavior and strength are concerned, the response at 
one material point and of a coupon under a uniform uniaxial load must be the same. 
However, that is not the case for the post-peak response due to damage localization. 
At the material point level, a stable post-peak strain softening must exist such that the 
area under the stress-strain curve is equal to the fracture energy per unit volume of the 
crack band. However, just like the experiments, a stable post-peak cannot be obtained 
in coupon level simulations. To ensure that the model dissipates the correct fracture 
energy, a material point calibration is required. 
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4-2.1. Calibration: mode 12 


As shown in the foregoing sections, mode 12 is strictly related to the in-plane nniaxial 
behavior of the material. Accordingly, the nniaxial tests on the [0°]8 specimens were 
nsed for the calibration of the 6 parameters reqnired to describe the inelastic behavior 
in tension and compression. The parameters atu and Ocu were 

calibrated in snch a way that: (1) the predicted strength matched that from the [0°]8 
conpon tests, for both tension and compression; (2) the predicted size effect in strength 


matched the test data (as shown in section 4.3.3). 

The single material point simulations were conducted with reduced integration on 
one shell element of size he = 2 mm. The calibrated parameters are presented in Table 
It should be noted that, even if the exponents ati 2 and aci 2 might have different values, 
an excellent agreement with experimental data was obtained setting 0*12 = aci 2 = 3/4. 

Following the material point calibration, a verification was performed with a coupon 
level simulation. For the uniaxial tensile test coupon, the dimensions were 100 mm x 
40 mm X 1.9 mm, and for the compressive test they were 100 mm x 40 mm x 1.9 mm, 
in accordance with the ASTM standards [5H1 [59] . The hnite element models for both 
coupons were meshed with four-noded shell elements of size 2 mm. 

The results for the coupon simulation are shown in Figures ([^,b). As expected, 
despite a stable post peak at the material point level, a sudden dynamic failure of 
the coupon is observed at the peak load. This can be seen as the vertical drop of 
load induced by snap-back instability. Both figures show excellent agreement with 
experiments under tension as well as compression. 


4 . 2 . 2 . Calibration: mode 4 

The behavior of the material under off-axis uniaxial loading is highly dominated by 
the matrix, in particular, by its in-plane shear response. In this conhguration, mode 4 is 
predominant so these tests can be used for its calibration. The mode 4 parameters 

and p were calibrated at material point level, to match accurately 
the pre-peak non-linearity and the peak load from the [45]8 coupon. The fracture energy 
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in this configuration was not known. Thus, a conservative value of about 50 N/mm was 
assumed and the parameters regulating the post-peak response were 

calibrated accordingly. 

In the absence of experimental data, the same set of values used in tension was used 
to calibrate the compressive response. It is worth mentioning, however, that the model 
does have the capability to describe a different behavior in tension and compression 
under off-axis loading, a useful feature to describe the shear response of nonsymmetric 
textile composites. 

Following the material point calibration, coupon level simulations were conducted 
for verification using a 2 mm mesh. The coupon dimensions, in accordance with the 
ASTM standards, were 100 mm x 40 mm x 1.9 mm. 

Figure ([^) shows the specimen load-displacement curves. It is seen that the shear 
response is very well matched by the model, including the pre-peak non-linearity caused 
by constrained matrix microcracking. An instability is observed at the peak load, 
manifested as a vertical drop of load to zero, which is indicative of a snap-back. The 
resulting calibrated values are shown in Table 

It is worth mentioning here again that the results in the 0° configuration were 
unaffected by the adjustment of the mode 4 parameters. This means that the model 
calibration can in general follow a clear sequence, avoiding the complication of an 
iterative procedure. 


4-2.3. Calibration: mode 3 and mode 5 


As pointed out in sections |2.4.4| and |2.4.6[ mode 3 and mode 5 are not expected to 
reach their inelastic boundaries as, in most practical cases, the material is considered 
to be almost in plane stress condition. Accordingly, in the absence of experiments 
to fully characterize these two modes accurately, their inelastic boundary can be just 
approximated. The parameters for mode 3 were estimated from [SSI EZ] such that the 
out-of-plane uniaxial elastic properties and strength match the one of a pure epoxy 
matrix. Mode 5 parameters were assumed to have the same values as for mode 4. The 
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calibrated parameters are shown in Table 

4 . 3 . Verification and validation 

4 . 3 . 1 . Uniaxial tests on Q.I. laminates 

Snbseqnently to the parameter calibration, predictions were carried ont to validate 
the model. The tests used for this purpose were uniaxial tests on coupons with a quasi¬ 
isotropic lay-up, [0°/45°/—45°/90°]s, for both tension and compression. The coupon 
dimensions were 100 mm x 40 mm x 1.9 mm in tension and 100 mm x 40 mm x 1.9 
mm in compression, in accordance with the ASTM standards. Finite element models 
of the coupons were built and meshed with layered shell elements of size 2 mm, 3 Gauss 
points being assigned to each layer. 

The comparison between experimental and predicted load-displacement curves is 
plotted in Figures ([^,b). In both cases, the agreement in the initial slope and strength 
is extremely good. For the tensile tests, the predicted strength is slightly on the lower 
side. This may be attributed to the fact that the mean strength measured from the 
0° conhguration tensile tests, used in calibration, had itself a signihcant scatter. The 
prediction of the tensile strength of the QI coupon rests heavily on this calibration. 
The reason is that the hrst layer to experience failure is indeed the layer oriented along 
the loading axis since the strain at failure for this orientation is low. The agreement 
for the compression case is excellent. 

It worth mentioning here again that these results were obtained without changing 
any of the parameters calibrated before. The excellent agreement with experimental 
data serves as hrm validation of the theory. 

4 . 3 . 2 . Biaxial failure envelopes 

In order to asses the capability of the model to capture the multi-axial behavior of 
the composite, the experimental multi-axial failure surfaces reported in [HTJ [69] were 
considered. It is worth mentioning here that the material system studied in these 
references is a plain weave composite, different from the twill 2x2 composite used to 
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calibrate the model. Accordingly, the comparison between model predictions and the 
experimental data has to be considered as qualitative only. 


Figure (lOr) shows the experimental data by Owen and Griffiths 0 on the failure 
under normal biaxial loading of a system composed by a polyester matrix reinforced 
by a glass woven fabric. In order to compare these results with the calibrated model, 
the failure stresses are normalized against the respective uniaxial strengths with 
i =2,3. As can be noted, the model shows an excellent agreement with experimental 
data, a remarkable result considering that the spectral stiffness microplane model does 
not require any additional parameter to seize the multi-axial behavior of the material. 
As a reference, the Tsai-Wu failure envelope [68], calibrated using only the uniaxial 
strengths in direction 2 and 3, is also plotted. As can be seen, the predictions from 
the microplane model and the Tsai-Wu criterion are very similar, with the microplane 
model being slightly more accurate especially in the II quadrant. 

The behavior of the material under tension-torsional loading is represented in Figure 
([IoId) which reports the experimental data by Fuji! et. al [Uni on tubes made of a 
polyester matrix reinforced by a glass plain woven fabric. It can be seen that the 
model compares satisfactorily with the experimental data even if it tends to slightly 
overpredict the resistance to failure. The hgure also shows the prediction provided by 
the Tsai-Wu criterion calibrated only with uniaxial data. 


4-3.3. Size effect tests 

The fracture energy that was indirectly measured from size effect tests was incor¬ 
porated in the model by calibrating the material point response. In-lieu of the crack 
band model, the calibration satishes the condition that the area under the stress-strain 
curve equals the fracture energy per unit volume (Figure [11)3). However, this does 
not guarantee that the model will be able to predict the right size effect in structural 
strength. This is because, for lab-scale structures of quasibrittle materials, the fracture 
process zone (FPZ) might not have reached its full size when the structure reached its 
peak load. Therefore, as was shown by Cusatis et al. ra. the peak load is governed 
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by not only the numerical value of the fracture energy, but also by the shape of the 
softening law. Thus, checking the predicted size effect is an important verihcation of 
the damage formulation. Material failure theories that are based purely on strength 
such as e.g. Tsai-Wu [HH] cannot make this prediction. This aspect will be considered 
in detail in a separate forthcoming article [HI]. Here, only the comparison between the 
predicted size effect and the available test data is shown. Further, the Size Effect Law 
(SEL) proposed by Bazant [21] htting the experimental data is reported as a reference. 

Finite element models for all three SENT coupons were built and uniaxial tensile 


simulations were performed to predict the peak load. As noted in section 4.1.1 , the 
fracture propagation for these specimens was modeled in the sense of the crack band 
model. The elements in the crack band had a height of 2 mm whereas the size in 
the direction of crack growth was 0.2 mm to capture correctly the stress prohle in the 


ligament. Figure (Hr) shows, for the coupons of different size, a comparison of the 
predicted and measured strengths, dehned as aj^ = P/tD with P =maximum load, 
t =thickness and D =witdh of the specimen. The calibrated material point stress- 


strain curve, instead, is depicted in Figure (Hr). It can be seen that, indeed, the 
model predicts the correct size effect. It should be noted that the specihc assumed 
form of the hber softening is instrumental in this prediction. Furthermore, using the 
size effect htting, the total fracture energy, Gp, estimated from the modeling predictions 
was 88.4 N/mm, which agrees very well with the value of the initial fracture energy, 
Gf, estimated from Size Effect Law 74 N/mm [61]). In the absence of experimental 
data, a similar value was assumed in compression. 

The capability of capturing intra-laminar size effect of textile composites is key for 
safe damage tolerance design across various engineering applications. This is dehnitively 
a feature of the spectral stiffness microplane model as the excellent agreement with 
experiments showed. 
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5. Application to crushing of composite tubes: results and discussion 


Thanks to their outstanding energy absorption capability, textile composites are 
often used as substitute of steel and aluminium to enhance crashworthiness of structural 
components. However, in order to fully exploit the potential of these materials, material 
models to be used as design tools are key. For this reason, the predictive capability of 
the spectral stiffness microplane model was challenged to predict the energy absorbed 
during the impact of composite crash cans. Figure dlTj^) shows the geometry of the 
structure under study consisting of a hat section tube and a reinforcing plate glued 
together by a toughened epoxy glue. Two different lay-up configurations were studied, 
namely: 


1. hat section tube: [0°]ii, plate: [0°]8; 

2. hat section tube and plate: [07907457-4570790707-4574579070°]; 


Following the test conditions described in section 3.4, the crush cans were modeled in 
Abaqus Explicit [62j using a mesh of triangular shell elements of 2 mm. All the degrees 


of freedom of the nodes at the bottom section were fixed while a predefined velocity 
field of 4.6 m/s was prescribed to the impacting mass of 74.4 kg consisting of rigid shell 
elements. The general contact algorithm provided by Abaqus Explicit |H2] was used 
while element deletion was adopted to avoid excessive element distortion during the 
simulation. The elements were deleted as soon as the dissipated energy reached 99% 
of fracture energy or if the magnitude of the maximum or minimum principal strains 
reached 0.45. 

The comparison between experimental and numerical results for case 1) and 2) is 
shown in Figures ([II^,b) in terms of reaction force on the plate versus time while Figure 
(fitj) shows the typical predicted fracturing pattern for case 1). As can be noted, a 


satisfactory agreement is found for both cases. As shown in Figure (Hr) for case 
1), the predicted plateau load, i.e. the reaction force on the plate once the crushing 
process is stabilized, compares very well with the experimental value: 33.7 kN vs 35.6 
kN. However, not the same accuracy is provided for the time required to stop the mass. 
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which is slightly overpredicted. This is not surprising since, as can be noted from Figure 


(Hr) the reaction force at the early stages of the simulation is underpredicted. This 
might be an artifact of element deletion in the first stages of contact between the plate 
and the crush can and it is not considered here related to the material model. This 
was confirmed by the fact that running the same simulation with different material 
models provided the same problem in the early stages of the simulation. It is worth 
mentioning here that, during the tests, only 1/4 of the crush can was crashed, being the 
initial velocity not high enough. With an higher initial velocity, the prediction of the 
model in terms of time to stop the impacting mass would have been far more accurate. 
This because the extent of the region in which progressive stable crushing was present 
(predicted by the model with an error of about 5%) would have been much larger. 

Similarly, Figure ([TT)3) shows the results for case 2). In this configuration, several 
layers of the material were loaded off-axis dissipating energy by extensive matrix mi¬ 
crocracking before the onset of the first intra-laminar cracks. This damage mechanisms 
is supposed to provide extra energy dissipation compared to case 1) and, as a matter 
of facts, the experimental plateau load was in this case higher: 41.0 kN. Thanks to the 
excellent capability of capturing matrix microcracking, as was shown in the previous 
sections, the spectral stiffness microplane model was able to predict this increased force 
with a value of 42.6 kN. This proves that the model can be used as design tool to find 
the optimum lay-up configuration to maximize energy dissipation. It should also be 
noted that, in this case, also the total time to stop the impacting mass was predicted 
with very high accuracy. This may be in part due to the fact that the initial under¬ 
prediction of the reaction force in the early stages of simulation was compensated by a 
slightly overpredicted plateau force. 

It should be highlighted here that the results shown in Figures ([Il^-d) represent 
a pure prediction based on the calibration and validation done through uniaxial tests 
as well as size effect tests only. No adjustment of any of the parameters of the model 
was done, making the results shown even more remarkable. This achievement was made 
possible by the accurate modeling of the main damaging mechanisms such as e.g. matrix 



microcraking and, most importantly, by the correct modeling of intra-laminar cracking 
through the introduction in the material model of a characteristic length scale. The 
latter mechanism of energy dissipation seems to be the most important for the range 
of impact velocities considered in the present work, in which the maximum strain rate 
is not supposed to be higher than 20 s“^. As reported from experiments, delamination 
seemed to be rather limited. 

6. Conclusions 

This contribution proposes a theoretical framework, called “Spectral Stiffness Mi¬ 
croplane Model”, to simulate the orthotropic stiffness, pre-peak nonlinearity, failure 
envelopes, and the post-peak softening and fracture of textile composites. Based on the 
results presented in this study, the following conclusions can be formulated: 

1. The microplane formulation with constitutive laws dehned on planes of several 
orientations within the mesostructure allows a sound and physically-based de¬ 
scription of the damage mechanisms occurring in textile composites whereas the 
spectral decomposition of the microplane strains and stresses provides a rigorous 
generalization of the formulation to anisotropy; 

2. Each eigenmode can be easily associated to the mechanical behavior of a particular 
constituent at the mesostructure and to a particular type of deformation. This 
makes it easy to take advantage of each eigenmode to dehne constitutive laws at 
microplane level targeting different damaging and fracturing mechanisms; 

3. Applied to a carbon twill 2x2 composite, the model showed excellent agreement 
with uniaxial tests in tension and compression for various lay-ups. In particular, 
the model captured the highly non-linear mechanical behavior under off-axis and 
shear loading, which is typical of textile composites and is characterized by dif¬ 
fuse subcritical microcracking of the polymer matrix. This feature is of utmost 
importance in all situations in which predicting energy absorption is key, such as 
in crashworthiness analyses; 
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4. Thanks to the use of microplanes representing the mesostructure of the material, 
the model captures the complex multi-axial behavior of textile composites without 
the need of any additional parameter; 

5. Different from strength-based criteria abundant in the literature, the formulation 
is endowed with a characteristic length through coupling with the crack band 
model. This ensures objective numerical analysis of softening damage and pre¬ 
vents spurious mesh sensitivity. Further, this is key to capture the intra-laminar 
size effect, a salient feature of composite structures. This aspect, too often over¬ 
looked in the literature on composites, is a determinant factor for damage toler¬ 
ance design of large composite structures; 

6. Compared with experimental results on the axial progressive crushing behavior of 
composite crush cans, the model provided excellent predictions for all the lay-up 
conhgurations under study. The reason is that the main damage mechanisms such 
as matrix microcracking and longitudinal intra-laminar cracking were captured 
correctly by the formulation; 

7. The model is computationally efficient and capable of analyzing the fracturing 
behavior of large composite structures, making it a valuable design tool for crash- 
worthiness applications. Further, it has sufficient generality to allow extensions 
to composites with more complex architectures, such as the hybrid woven com¬ 
posites, and 2D or 3D woven or braided composites; 

8. The presented model is currently implemented in the commercial code MARS 
with the name “Woven Composite Lamina”. 
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Description 

Symbol (units) 

Measured value 

In-plane modulus 

E=E2=E3 (GPa) 

53.5 

Out-of-plane modulus 

E' = El (GPa) 

11.0^ 

In-plane shear modulus 

G = G23 (GPa) 

4.5 

Out-of-plane shear modulus 

G' = Gi 2 = Gi 3 (GPa) 

3.6^ 

In-plane Poisson ratio 

U = 1/23 = ^232 (-) 

0.055 

Out-of-plane Poisson ratio 

12 ' = Z/31 =1/21 (-) 

0.4^ 


^properties estimated from [561 EZ] • 


Table 1: Experimental elastic properties of carbon twill 2x2 composite 


Mode 


M) 




1 

2 

3 

4 

5 


a [cos^ ^ + (x cos^ If + sin^ fj sin^ 6] 

/3 (sin^ 9 — sin^ f sin^ 6) 

7 [cos^ 0 + (C f) sin^ 6] 

£4 sin f sin 29 


(x — 1)^ cos^ f sin^ 0(3 — cos 2f + 2 cos^ f cos 29) 


/3^ sin^ 9 sin^ 2f + cos^ 9 (cos 2f — 3)^ 

5 (C ~ 1)^ cos^ f sin^ 6(3 — cos 2f + 2 cos^ f cos 29) 
y^£|~(cos^~^'cos^~0~+Zos^~20'sin^~(^ 

£5 cos f sin 29 + eg sin^ 9 sin 2f (eg cos 9 sin f — eg cos 2f sin 9)^ + cos^ f (£5 cos 29 + ea sin f sin 29)^ 


where a = (£2 + £3 +eix) /(2 + X^)> /3 = 1/2 (£3 - £ 2 ), 7 = +£3 +eiC) /(2 + C^) 

Table 2: Microplane eigenstrains: normal and tangential components 
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Mode Description 

Symbol (units) 

Calibrated value 

mode 1 elastic eigenvalue 

mode 2 elastic eigenvalue 

microplane peak stress in tension 

parameter governing post-peak softening in tension 

12 

parameter governing post-peak softening in tension 

microplane peak stress in compression 

parameter governing post-peak softening in compression 

parameter governing post-peak softening in compression 

A(b (GPa) 

A(2) (GPa) 
(MPa) 

(-) 

a.i2t (-) 

(MPa) 

(-) 

a,ci2 (-) 

61.85'^ 

50.71'^ 

400 

30.60 X 10-3 

0.75 

405 

30.60 X 10-3 

0.75 

mode 4 elastic eigenvalue 

microplane stress in tension at start of non-linear boundary 

exponent governing pre-peak non-linearity in tension and compression 

strain at starting of post-peak softening in tension 

4 

parameter governing post-peak softening in tension 

microplane stress in compression at start of non-linear boundary 

strain at starting of post-peak softening in compression 

parameter governing post-peak softening in compression 

A(4) (GPa) 
(MPa) 

P (-) 

(-) 

ki? (-) 

4"*^ (MPa) 

fcit' (-) 

(-) 

8 .10^ 

45 

0.3 

124.6 X 10-3 

120.15 X 10-3 

45 

124.6 X 10-3 

120.15 X 10-3 

Calculated by means of Eqs. (|l2^-d). 



Table 3: Material model parameters, calibrated by means of uniaxial and 

size effect tests, 

describing the 

in-plane behavior of the composite. 
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Mode 

Description 

Symbol (units) 

Calibrated value 


mode 3 elastic eigenvalue 

A(3) (GPa) 

10.82'^ 


microplane peak stress in tension 

(MPa) 

90 


strain at starting of post-peak softening in tension 

fci?’ (-) 

4.0 X 10-3 

3 

parameter governing post-peak softening in tension 

(-) 

20 X 10-3 


microplane peak stress in compression 

(MPa) 

90 


strain at starting of post-peak softening in compression 

kif (-) 

4.0 X 10-3 


parameter governing post-peak softening in compression 

(-) 

20 X 10-3 


mode 5 elastic eigenvalue 

A(5) (GPa) 

7.20^ 


microplane stress in tension at start of non-linear boundary 

(MPa) 

45 


exponent governing pre-peak non-linearity in tension and compression 

P (-) 

0.3 


strain at starting of post-peak softening in tension 

(-) 

124.6 X 10-3 

5 

parameter governing post-peak softening in tension 

(-) 

120.15 X 10-3 


microplane stress in compression at start of non-linear boundary 

(MPa) 

45 


strain at starting of post-peak softening in compression 

(-) 

124.6 X 10-3 


parameter governing post-peak softening in compression 

4” (-) 

120.15 X 10-3 


^ Calculated by means of Eqs. ( jl2p -d). 

Table 4: Material model parameters related to the out-of-plane behavior of the composite. The parameters 
for mode 3 were estimated from [561157| . Mode 5 parameters were assumed to have the same values as for 
mode 4. 
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a) b) 


Figure 1: Schematic representation of a) the Representative Unit Cell of a 2x2 twill composite with its 
local coordinate system and the microplanes used to define the constitutive laws of the material; b) local 
spherical coordinate system. 
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Figure 2: Distribution of the normal microplane strain related to a) mode 1, b) mode 2 and c) model2 
in the presence of a uniaxial macroscopic strain in plane 2-3. 
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Mode 3 



Figure 3: Distribution of the normal microplane strain related to mode 3 in the presence of a uniaxial 
macroscopic strain ei. 
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Figure 4: Distribution of the normal microplane strain related to mode 4 in the presence of a macroscopic 
in-plane shear strain 64 . 
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Mode 5 



Figure 5: 


Distribution of the normal microplane strain related to mode 5 


in the presence of macroscopic 


out-of-plane shear strains and £g. 
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Figure 6: Spectral stiffness microplane constitutive laws, a) Typical stress boundary for mode 12 with 
exponential decay; b) Typical bilinear boundary for mode 3; c) mode 4 and 5 inelastic stress boundaries; d) 
typical unloading-reloading rule for a /-th generic mode. 
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Load [kN] 
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Displacement [mm] 

Figure 7: Comparison of the measured and calibrated load displacement curve for a) a [Ojs coupon, under 
tension; b) a [Ojg coupon, under compression and c) a [45]8 coupon, under tension. 



Displacement [mm] 


Figure 8: Comparison of the measured and predicted load displacement curve for a quasi-isotropic coupon 
[0°/45°/—45°/90°]2 s under (a) tension, (b) compression. 
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Normalized strength, log ctj 


a) b) 
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Figure 9: a) Comparison of the measured and predicted size effect in geometrically scaled Single Edge 
Notch Tension [Ojs coupons; b) typical stress-strain curve in pure tension. The total fracture energy Gp is 
calibrated adjusting the post-peak softening response of the material. 
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a) b) 



Normalized strength 


Figure 10: Comparison between measured and predicted biaxial behavior, a) Biaxial loading; b) tension- 
torsion loading. 
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Figure 11: Crashing of composite tubes, a) comparison between model prediction and experimental results 
(hat section: plate: [0°]8); b) comparison between model prediction and experimental results (hat 

section and plate: [0°/90°/45°/—45°/0°/90°/0°/—45°/45°/90°/0°]); c) geometric specification of the 


crash can cross-section (dimensions in mm); d) typical predicted fracturing pattern. 
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